Winter climate is a key factor affecting population dynamics in high-arctic ungulates, with many studies showing a strong negative correlation of winter precipitation to fluctuations in population growth rate. Terrestrial ice crust or ground-ice can also have a catastrophic impact on populations, although it is rarely quantified. We assess the impact of winter climate on the population dynamics of an isolated herd of Svalbard reindeer near Ny-Ålesund with a retrospective analysis of past winter snowpack. We model landscape-scale snowpack and ground-ice thickness using basic temperature and precipitation data in a simple degree-day model containing four adjustable parameters. Parameter values are found that lead to model snow and ground-ice thicknesses which correlate well with three different model targets: reindeer population growth rates; April snow accumulation measurements on two local glaciers; and a limited number of ground-icing observations. We explain a significant percentage (80%) of the variance in the observed reindeer population growth rate using just the modeled mean winter ground-ice thickness in a simple regression. Adding other explanatory parameters, such as modeled mean winter snowpack thickness or previous years' population size does not much improve the regression relation.
Introduction
Population dynamical studies often show that observed fluctuations in large herbivore population size arise from both density dependence and environmental stochasticity (Sæther, 1997; Grenfell et al., 1998; Gaillard et al., 1998, 2000; Post and Stenseth 1999; Coulson et al., 2001). In both temperate and arctic areas, winter climate has often been found to influence large herbivore population dynamics, commonly with a negative relationship between amount of snow and population parameters such as growth rate, survival or fecundity (Fancy and White, 1985; Adamczewski et al., 1988; Sæther, 1997; Gaillard et al., 1998; Post and Stenseth, 1999; Post and Forchhammer, 2003).
This pattern also includes Svalbard reindeer (Rangifer tarandus platyrhynchus) where both density dependence and a strong negative relationship between population growth rate and winter precipitation amount have been identified (Aanes et al., 2000, 2003; Solberg et al., 2001). Snow both reduces access to available food and increases energetic costs. This may reduce female fecundity and increase the mortality, especially for the youngest individuals (see Sæther, 1997; Gaillard et al., 1998 for reviews).
A less-studied climatic phenomenon is the formation of terrestrial ice-crust or ground-ice. It has long been appreciated that ground-icing events have a strong impact on population dynamics since they can completely block access to the food eaten by large herbivores. Formation of ground-ice, leading to a decrease in population size, or even to population crashes, has been reported both on Svalbard (Lønø, 1959; Alendal and Byrkjedal, 1974; Larsen, 1976; Reimers, 1977, 1982, 1983; Tyler, 1987; Øritsland, 1998), as well as in other areas (e.g., Parker et al., 1975; Forchhammer and Boertmann, 1993). Few of these studies have actually tried to systematically quantify the impact of icing on the dynamics; most studies report single observations of population decreases after severe winters. One exception is from Greenland, where Forchhammer and Boertman (1993) calculated an index to reflect the formation of ground-ice and evaluated how it affected muskox populations. Another example comes from Svalbard where Solberg et al. (2001) found a negative relationship between number of winter warm days (>0°C) with precipitation and Svalbard reindeer population growth rate. Finally, Putkonen and Roe (2003) look explicitly at the effect of rain-on-snow events, the main contributing factor leading to ground-icing. However, while they discuss qualitatively the effect this may have on biota, they did not include any of the available ecological data, on the local reindeer population, for example.
In the present paper we use information from a well-studied Svalbard reindeer population on the Brøgger peninsula area, near Ny-Ålesund (Fig. 1). A previous study found no evidence of strong density dependence operating in this population (Aanes et al., 2000). The strongest correlation was between the population growth rate (G) (by custom referred to as Rt in the biological literature), given by
a measure of the change in population size N between years t and t − 1, and the summed winter (Oct.–Apr.) precipitation (Aanes et al., 2000). A dramatic effect of winter precipitation on this reindeer population is best exemplified by the ground-icing event that occurred on the Brøgger peninsula in autumn 1993. During an 8-d period in November and December, record amounts of rain fell, while air temperatures fluctuated around 0°C. Prior to the rainfall, there had been a prolonged cold period, with no snow, and ground temperatures were well below freezing. The net result was a thick ice layer (up to 10–30 cm) deposited over large parts of the lowland plain, and a substantial reduction in the reindeer population the following spring. While the population reduction was not wholly due to mortality (there was migration over the glaciers and frozen bays to previously uncolonized areas), a substantial number of reindeer died in the winter of 1993/94 due to the icing.
Given the apparent effect of winter precipitation on population change, we would expect that the relation between population growth rate G and winter precipitation is improved when one uses actual snowpack thickness data. Furthermore, there may have been other, less dramatic episodes of icing in the past, which were not noticed or deemed noteworthy, but which nonetheless restricted reindeer's access to food. Even the most significant icing event, which is only reported orally (Aanes et al., 2002), was not measured systematically in the field, and until recently, there have been no regular observations of the snowpack in the area. Our only recourse, then, is to model the snowpack with available meteorological measurements. Using a simple snowpack model driven by twice-daily precipitation and temperature measurements, we retrodict snowpack thickness and the formation of ground-ice, and evaluate its importance for Svalbard reindeer dynamics in the area.
Field Site and Data
FIELD SITE
The Brøgger peninsula, located on the northwestern coast of Svalbard, is a rugged, glaciated mountain ridge, surrounded by an extensive lowland plain (Fig. 1). Total area is about 220 km2, of which 25% is glaciated and 50% contains 5% or greater vegetation coverage. Vegetation comprises small vascular plants, grasses, and cryptogams (Øritsland et al., 1980; Brattbakk, 1986). The climate is characterized by low mean annual temperature (−6°C), frequent winter warm temperature excursions, and low mean annual precipitation (400 mm).
REINDEER
There were no Svalbard reindeer on the peninsula from the beginning of the 1900s, when the population was eradicated by hunting, until 1978, when 15 reindeer were introduced as part of a long-term ecology experiment. Open water and glacier fronts provide natural barriers to impede migration of the introduced population, under normal conditions (Fig. 1). From 1979 to the present, the total number of reindeer on Brøgger peninsula has been counted each spring (Aanes et al., 2000); these measurements are the basis for the annual population used in this paper (Fig. 2).
METEOROLOGICAL MEASUREMENTS
Meteorological data from Ny-Ålesund for the period 1 January 1969–31 December 2002 are taken from routine measurements made by the Norwegian Meteorological Institute. Input data are twice-daily (0700 and 1900) measurements of standard 2-m screen air temperature and summed 12-h precipitation from a heated rain gauge. The twice-daily temperatures are averaged to estimate the effective temperature (T) of the preceding 12-h period during which time an amount of precipitation (P) collected in the rain gauge.
Precipitation data are not adjusted to account for wind speed, which reduces the amount of precipitation caught by the gauge, particularly for snow (Førland and Hanssen-Bauer, 2000), nor do we apply any elevational correction factor (Førland et al., 1997), since we have found that the model results are insensitive to either precipitation correction factors.
SNOW AND ICE MEASUREMENTS
Very little is recorded about ground-icing, either on Svalbard in general or in Ny-Ålesund in particular. A summary of all known or measured icings on the Brøgger peninsula is given in Table 1. There are only sporadic measurements made of snow depth at the meteorological site in Ny-Ålesund proper. Long-term annual snow-depth measurements are available in the area from the winter glacier mass-balance records (e.g., LeFauconnier et al., 1999) of the nearby glaciers Midre Lovénbreen and Austre Brøggerbreen (Fig. 1). The mass-balance data used here (J. Kohler, unpublished data) cover the period 1970–2002, and comprise averages of snow-thickness measurements made in late April–early May at the lowest elevations (50–150 m). While restricted spatially to the two glaciers, the data (Fig. 3) are good proxies for spring snowfall over the entire Brøgger peninsula, as demonstrated by parallel measurements made over the period 2000–2003.
Model
MODEL DESCRIPTION
We use a simple degree-day model (e.g., Rango and Martinec, 1995) to simulate Ny-Ålesund snowpack and icing back in time. The model is essentially a nonlinear filter operating on time-series of temperature (T) and precipitation (P). The filter highlights the thresholds and interactions between T and P that provide the fundamental controls on whether or not ground-icings can form. Our model is not a “physical” snowpack model, such as SNTHERM (Jordan, 1991) or SNOWPACK (Bartelt and Lehning, 2002), which can simulate the detailed layer structure within a time-varying snowpack. Physical models are inappropriate for the problem at hand for several reasons. First, they include a large number of relevant but difficult to model processes (e.g., snowpack densification, snow crystal metamorphosis, and water transport) which require in situ measurements to parameterize or verify. Second, the input and verification data are difficult to obtain as continuous time-series since they generally require precise calibration and a high degree of homogenization (e.g., radiation, albedo, wind profiles). Third, and more importantly, the results from physical models still need to be extrapolated to the larger landscape or regional scale from the location(s) in which the required data are available. Since in most cases the extrapolation is as poorly constrained for physical models as it is for simpler models, it ultimately comprises the largest source of mismatch between model and reality (Rango and Martinec, 1995).
The model (Fig. 4) comprises two boxes: a snowpack box and a ground-ice box. Precipitation occurs either as rain or as snow; snow is added to the snowpack box, while rain or water formed from melting in the snowpack box is assumed to pass directly into the ground-ice box. Melting of snow or ground-ice follows the degree-day concept, that is, the amount of melt is indexed to the number of degrees per day above a certain threshold (Rango and Martinec, 1995).
Four parameters control the behavior of the model: two threshold temperatures, TSR and TDD0, the first of which determines the precipitation type and the second of which determines the temperature above which snow can melt, and two degree-day coefficients, DDS and DDI, which control the amount of snow and ice melt.
Input data are air temperature T(t) and precipitation P(t) as a function of time t. At each time-step ti, when the air temperature T(ti) is less than a threshold temperature TSR, precipitation occurs as snow, and the snowpack thickness S(ti) increases by the precipitation amount P(ti). When T(ti) is greater than or equal to TSR, precipitation occurs as rain. Rain is assumed to pass through the snowpack box to fill the ground-ice box by the precipitation amount P(ti), if the ground is frozen (we assume this to be the case between the months of October and May; see below for further discussion).
When T(ti) is greater than a second threshold temperature TDD0, the snowpack thickness S(ti) is reduced at a rate DDS [T(ti) − TDD0]. Melted snow is treated like rain, as described above, moving to the ground-ice box when the ground is frozen. If the snowpack is absent then the ground-ice melts at a second rate DDI [T(ti) − TDD0].
Melting ice, or any melt or rain falling between June and September, passes out of the model. Snow or ice remaining in the model at 1 August is removed.
MODEL ASSUMPTIONS
Snowpack density is not considered in the model; water formed in the snowpack simply passes directly to the bottom to form ground-ice, neglecting the possibility that ice or high density layers form within the snowpack. For the problem at hand, modeling the effect of winter snowpacks on reindeer, this is not a completely unreasonable simplification since the formation of ice layers anywhere within the snowpack provides a similar, albeit more easily breached, physical barrier to the food source.
We ignore the complication of modeling the heat budget of the ground and the snow layer (e.g., Putkonen and Roe, 2003), since this introduces new assumptions and parameters concerning land-surface, soil, and permafrost properties. Instead we simply assume that the ground is frozen between the months of October and May. Records of ground temperature from the Ny-Ålesund area show that this assumption is not unreasonable (Putkonen, 1998; Roth and Boike, 2001; Boike et al., 2003; S. Gerland, unpublished data), although Putkonen and Roe (2003) show that the latent heat released in freezing during one major mid-winter icing event on the Brøgger peninsula raised the ground temperature to the freezing point. Nevertheless, this effect is short term and furthermore may not reflect the temperature at the icing surface.
Finally, we remove all remaining snow and ice on 1 August, if necessary. This is functionally equivalent to allowing higher summer ice melt rates, or in other words, temporally varying DDI. Degree-day factors do in fact vary seasonally (e.g., Rango and Martinec, 1995), but accounting for this in the model adds to model complexity. In any case, our application of the model in this paper is directed to winter effects. At this point we have too few data to constrain the behavior of the model in summer, apart from our observations of an ice-free landscape in all summers we have been there.
PARAMETER CHOICE
Landscape-scale models that use the degree-day approach (e.g., Rango and Martinec, 1995; Lindström et al., 1997) typically need to be calibrated, that is, model parameters must be adjusted to get the best fit of observed to modeled data. Parameter values reported in the literature, particularly degree-day coefficients, can vary widely (e.g., Braithwaite, 1995; Rango and Martinec, 1995; Hock, 1999, 2003; LeFebre et al., 2002). Therefore, we cannot simply assume that the parameters take on any reported values.
We explore parameter space to see how different values affect the model results. We assess model goodness-of-fit using the following targets: (1) the reindeer log population growth rate G (Fig. 1), which we assume to be an ground-icing proxy (2) the winter mass-balance measurements (Fig. 3), and (3) the few observations of icing (Table 1).
It may seem imprudent to use reindeer population growth rate as a ground-icing model target, without considering any other internal population factors, but given that previous work has demonstrated a clear connection between winter precipitation and population change (Aanes et al., 2000, 2003), and that two of the major population downturns occurred in years for which we do know there was icing (Table 1), it is not entirely unwarranted to use changes in the reindeer population as a ground-icing indicator. Furthermore, there are simply too few measurements of ground-icings (Table 1) to use these data exclusively as a model target. Finally, it is not possible to use the mass-balance data as a ground-icing model target since these data are snow depths, that is, the thickness obtained from probing the snowpack until something hard is reached. Whether this is ground or ground-ice is not recorded, so that in this case the comparison of modeled to observed data would be completely insensitive to the choice of parameters that suppress or promote model ground-icing.
MODEL CALIBRATION
The output data from the model are time-series of twice-daily snowpack thickness S, and ground-ice thickness I. We assume an expected range of values for the model's four adjustable parameters and, moving systematically at reasonable step-sizes between these minimum and maximum values (Table 2), obtain trial S* and I* time-series for each combination of parameters. The trial model results are then compared quantitatively to targets via three regressions.
REGRESSION 1
The first (hereafter referred to as R1) is a multiple linear regression that assumes that reindeer log population growth can be modeled simply by
Here Sw* is the average winter (Oct.–Apr.) snowpack thickness, Iw* is the average winter (Oct.–Apr.) ground-ice thickness, and both quantities are derived from model S* and I* trial time-series for the years 1979–2002. Hatted variables denote a regression result. The coefficients ci in equation (2) are determined from the regression, whose goodness-of-fit is evaluated by the familiar regression coefficient r2 (Davis, 2002).
Thus for each combination of the four adjustable model parameters, we obtain a value of r2 which indicates how well the Sw* and Iw* values obtained with those parameters did at explaining the observed reindeer population growth rate for the years 1979–2002.
Squaring the Īw term in equation (2) is a simple, low parameter cost way of expressing the nonlinear threshold effect of ground-ice on reindeer population dynamics. Ground-ice should have a threshold effect because thin ice layers provide little resistance to the sharp hooves of Svalbard reindeer. Above a certain thickness, however, there is an increasing energy cost to break through. Furthermore, the thicker the ground-ice layer is at an arbitrary point, the greater is the portion of the landscape covered by ground-ice.
REGRESSION 2
The second regression (R2) compares model 1 May snowpack thickness SM1* to bw100, the winter mass-balance data for the period 1970–2002 (Fig. 3) via
where c4 is a coefficient determined in the regression. In R2 we force the regression to pass through the origin. We also assume that winter mass balance was measured each year on 1 May. The actual dates vary by as much as 2 wk around this date, but are essentially unknown for much of the measurement period. This would contribute only small differences between model and observed data, however, since snowpack thickness is highly correlated over such short timescales during this period.
REGRESSION 3
The third regression (R3) compares IM, the four observations of icing thickness (Table 1), to the simulated May 1 ground-ice thickness IM1* via
where c5 is a coefficient determined in the regression. Data from 1995/96 is a qualitative observation, so we assume a value for that year (Table 1).
Regression Results
For each regression R1-R3, goodness-of-fit coefficients (r2) are computed for every combination of the four model parameters TDDO TSR DDS, DDI. Figure 5 shows plots of the r2 surfaces as a function of TSR and DDS, organized by regression number (R1-R3) and different values of TDDO; multiple surfaces within one plot correspond to different values for DDI. Figure 6 shows histograms of values for the occurrence of the four model parameters for >95th percentile results in each individual regression R1-R3, about 300 model results for each regression.
BEST VALUES FOR DDS AND DDI
In the R1 and R3 plots (Fig. 5), surfaces with lower regression coefficients generally correspond to higher values of DDI. In the R2 plots there is only one surface since the simulations lead to identical results regardless of the choice of DDI, not surprising given that R2 is based only on observed snowpack thickness. There is little difference, however, in R1 and R3 for different values of DDI at TDDO ≥ 2, the area of parameter space in which all three regressions have maximum r2.
The best regression results for DDS (Fig. 6) tend to lie toward values (≥10 mm °C−1) higher than those reported in Hock (2003) for a variety of mountainous settings, although Lefebre et al. (2002) find values as high as 14 mm °C−1.
In contrast, the best fits in R1 are associated with lower values of DDI (Fig. 6). Values of DDI are typically about 50% greater than DDS (e.g., Hock, 2003); however, although lower values of DDI appear to be preferred by R1, higher values lead to nearly equally good fits, as judged from the relatively close-spaced surfaces for TDDO > 1 (Fig. 5). Furthermore, the flatness of surfaces in areas with the highest r2 values (Fig. 5) shows that there are no truly preferred values for DDS.
BEST VALUES FOR TDDO
R1 and R2 both show the highest values for TDDO > 0°C. R3 has two maxima, one at low values of TDDO (−3 to −2°C) and the other at higher values (Fig. 5). Low values of TDDO result in more modeled ground-ice due to an increased incidence of melted snow; however, this does not necessarily reflect reality given that R3 is based on only four points, two of which are poorly constrained. Furthermore, the high r2 values at low values of TDDO are associated with unrealistically low values of DDI. Increasing DDI leads to a substantial reduction in the goodness of fit.
Accordingly, we assume that despite the first maximum in R3, all of the regression results consistently indicate only positive value for TDDO. For R1 and R2, best fits are obtained with TDDO at around 3°C, while in R3, best fits occur at 2 to 5°C (Fig. 6). In the literature (e.g., Rango and Martinec, 1995) TDDO is generally taken to be 0°C, so that starting melt only for temperatures higher than 3°C or more may seem unduly high. Part of the positive offset (ca. 0.5°C) could be explained by the temperature lapse-rate, since the model driving data are taken at near sea level while the observed data used in the regressions are effectively from higher elevations. Another possibility is that while snowmelt does in fact occur when temperatures are around 0°C, this water typically refreezes in the snowpack. Snowmelt only goes to bottom of snowpack to form ground-icing when higher temperatures are reached. A final possibility concerns seasonality in degree-day factors (Rango and Martinec, 1995). The combination of relatively high values of TDDO and constant degree-day factors serves to reduce winter melting and increase summer melting, duplicating the usual pattern for temporally-varying DDS and DDI (Rango and Martinec, 1995).
BEST VALUES FOR TSR
In R1 the effect of TSR is fairly weak (Fig. 5), although there is a slight preference for TSR = 0°C (Fig. 6), while in R3, there is a clear maximum for TSR = 0°C (Figs. 4, 5). The best fits in R2, on the other hand, occur predominantly at TSR = 2°C (Fig. 6). In either case, however, the obtained values are close to that indicated by the relation between observed precipitation type and 2-m temperature, whose transition occurs over the range 0–2°C (Fig. 7) as judged from the available data at Ny-Ålesund of observed precipitation type grouped by temperature at the observation time. The difference probably relates to the quantities being regressed; R1-R3 primarily respond to modeled ground-ice, which give better results with more rain (lower TSR), while R2 gives best results when there is more snow (higher TSR).
Best-Guess Model Snowpack
To summarize, model calibration shows that the good fits in R1-R3 occur over a wide range of parameter values, and that, with the possible exception of TDDO, these values are at least in the vicinity of expected or previously reported values. Accordingly, we now model a time-series of the “best-guess” snowpack thickness SBG and ground-icing thickness IBG (Fig. 8) using parameter values (TDDO = 3°C, TSR = 1°C, DDS = 8 mm °C−1, DDI = 12 mm °C−1) which could be taken from the literature, but which are guided to some extent by the results of R1-R3. The model is relatively insensitive to the ultimate choice of adjustable parameters, however, as demonstrated by comparing SBG and IBG to average thicknesses obtained from the 300 best simulations, the latter having been obtained by ranking the sum of the normalized r2 score in R1-R3. Qualitatively, there appears to be little difference between the two simulations (Fig. 8). Comparing observed G to the mean winter (Oct.–Apr.) ground-ice thickness ĪwBG, as well as those obtained from the 300 best simulations, shows that the two have consistent values for most years (Fig. 9). This is likewise the case for a comparison of observed mass balance to the mean 1 May snowpack thickness SM1BG computed from the best-guess and the 300 best simulations (Fig. 10).
Figure 10 also shows that modeled snowpack thicknesses are roughly half of those actually observed. This may be due to a number of factors, including spatial variability, elevation difference (Førland et al. [1997] report a ca. 20–25% increase in precipitation per 100 m elevation), and wind (which not only reduces precipitation catch in the gauge but are often stronger in Ny-Ålesund than on the glaciers).
Snow, Ice, and Reindeer Dynamics
We now apply two annual ground-ice and snow parameters [ĪwBG]2 and SM1BG derived above. We regress the observed reindeer population growth rate G against different combinations of these two modeled parameters as well as two observed parameters, summed winter (Oct.–Mar.) precipitation ∑ Pw, and log of the previous years' population ln N(t − 1). All regressions include a constant term. The parameters used in each regression and the resultant statistics are given in Table 3. Note that we do not combine regressions involving ∑ Pw and the modeled parameters since they are strongly correlated. Regressions are performed for all data from 1978/79 to 2001/02 (Table 3a), for all data but removing the extreme ground-icing year 1993/94 (Table 3b), and for all data but removing both 1993/94 and 1995/96 (Table 3c). The latter regressions are to test the robustness of the results when removing the extreme events.
As expected, there is a negative relationship between G and all winter climate parameters. And as reported previously (Aanes et al., 2000), the regression of G and ∑ Pw explains a significant amount of variance. But when considering all of the data, the square of the ground-icing parameter [ĪwBG]2 alone (regression 3, Table 3a) is by far the best single parameter for explain fluctuations in reindeer population growth rate, explaining more than twice the variance as explained by ∑ Pw alone (regression 1, Table 3a). Adding additional parameters in the regressions (e.g., regressions 5 and 6, Table 3a) does not improve greatly the amount of variance explained.
A time-plot of G obtained from regression of the ground-icing parameter shows that it captures the dynamics reasonably well (Fig. 11). There is, however, a temporal trend to the residuals (r2 = 0.28), with an increasing tendency toward over-prediction of the population growth rate with time. Although this temporal trend is dependent entirely on the first two and last three points, a likely explanation could lie in the grazing history of Brøgger peninsula. The reindeer were introduced to an area with superabundant resources, since there had been no grazing for roughly 100 yr. After the increase and crash in the reindeer population, studies have shown that the vegetation has been both reduced in amount and altered in composition (e.g., Cooper and Wookey, 2001; Henriksen et al., 2003). Hence, there has been an increasing imbalance between available food and reindeer density, and extreme climatic events such as ground-icings could be having an increasingly deleterious effect on the reindeer population.
Including the snowpack parameter SM1BG (regression 5, Table 3a) does not alter the basic result obtained using the ground-ice parameter alone; there is a slight improvement in the goodness-of-fit, but the temporal trend in the residuals (not shown) remains. However, while including the previous year reindeer density does not increase the goodness-of-fit by much (regression 5, Table 3a), it does reduce the temporal trend in the residuals (Fig. 11). While the temporal trend is highly dependent on a few data-points at the beginning and end of the time-series, we would argue nonetheless that the existence of a trend shows that the importance of population density has increased with time. This is to be expected in introduced populations following their initial growth phase with abundant food supply. Strong density dependence likely resulting from food competition has been found in two other Svalbard reindeer populations (Aanes et al., 2003), and we expect the Brøgger population to more or less mirror this dynamical pattern in the future, that is, with density dependence and environmental stochasticity operating in concert to influence the population dynamics.
There are three final points to consider. The first is that since ground-icing and G appear to resemble one another mostly on the basis of just two data points, it is natural to suspect that these extreme events alone are responsible for the high values of r2. Removing one or both years with heavy icing from the regressions (Table 3b, 3c) does result in a drop in the amount of variance explained by the ground-icing parameter. In fact, for the case in which we remove both extreme years, any of the parameters used alone (except snowpack) explain a similar amount of variance in G. Now the full regression (regression 6, Table 3c) explains more variance than the ground-icing parameter alone, implying that there is some underlying density-dependence in the population dynamics (see above). In addition, while r2 values decline in all regressions when removing the extreme years, the regression constants (not shown) obtained do not vary by much, suggesting that the regressions relations between the various parameters and G are relatively robust and not overly guided by outlier values.
The second point is that the ground-icing parameter [ĪwBG]2 was developed, in part, using precisely a regression involving the reindeer population (regression 5, Table 3a is identical to equation 2), so it may seem like circular reasoning to commend the snowpack model on its performance in predicting reindeer population growth. We have demonstrated, however, that the snowpack model results are insensitive to choice of model adjustable parameters. Furthermore we obtain similar results selecting the adjustable model parameters based on the R2 and R3 results alone.
Finally, as mentioned previously, squaring the ground-ice term ĪwBG is a simple, low parameter cost way of expressing the nonlinear threshold effect of ground-ice on reindeer population dynamics. The exponent is arbitrarily chosen to give a parabolic relation; improved fits to the data can be obtained by adjusting the exponent (in this case the best fit exponent for the regression is 2.6). While the data series is too short to explore the nature of such threshold patterns in detail, we can reach the same fundamental conclusions in this section with just a linear regression. In the latter case the regression coefficients are lower (e.g. Fig. 9) when using all of the data (i.e. Table 3a), but are essentially the same for the censored regressions.
Conclusion
Previous studies on Svalbard reindeer dynamics have shown that a relatively large proportion of the variance in population growth rate fluctuations can be explained from the amount of winter precipitation (Aanes et al., 2000, 2003). In this paper we refine our understanding of the importance of winter climate by decomposing winter precipitation amounts into snow and ground-ice components in the study area, using a simple modeling approach, and show that simulated ground-ice thickness is by far the best single parameter to explain the variance in reindeer population growth rate.
The ideal on study arctic population would include long time-series of in situ data on snowdepth, quality of snow (e.g., ice layers in the snow) and ground-ice thickness; such measurements are rarely, if ever, available to field ecologists, at least not over longer time spans. The simple modeling approach outlined in this paper provides a viable alternative given that standard measurements of precipitation and temperature usually available from meteorological stations close to a given study area. Our results show that the simple model is adequate to explain snowpack effects when evaluating mechanisms underlying often-observed fluctuations in large herbivores, with special emphasis to arctic areas.
The simplicity of the model allows rapid exploration of adjustable parameters, which may be necessary when extrapolating meteorological measurements made at one site to an entire landscape or region. This points to another direction forward, which is to consider the spatial or elevational dependence of snowpack characteristics. For example, we have observed reindeer foraging at high altitudes in the steep mountainous areas of Brøgger peninsula during winters with deep snow, or with ground-ice. It also remains to investigate what is available in these higher foraging grounds. The vegetation in the area is commonly mapped below the 200-m contour line (Brattbakk, 1986). If available forage above this height mainly consists of lichens it is plausible to believe that there are no large quantities, and the regrowth of lichens are extremely slow. Hence, such foraging refuges may be overgrazed rapidly if repeated blocking of food on the lowland plain happens.
Finally, another step forward in implementing the snowpack model is to drive the model with future climate scenarios. For example, a simple population dynamics relation of the type implied by the regressions summarized in Table 3 would then allow prediction of the future for Svalbard reindeer. As pointed out by Putkonen and Roe (2003), the frequency of climatic conditions which result in winter ground-ice may increase substantially in the future as climate shifts toward warmer and wetter winters, suggesting that more ground-icings are likely and that the future may be challenging for Svalbard reindeer.
Acknowledgments
We are indebted to Nils Are Øritsland who conceived and implemented the introduced reindeer experiment at Brøgger peninsula. We would further like to thank Erling J. Solberg and Regine Hock for comments which substantially improved the organization and readability of the paper. We are grateful to the Norwegian Meteorological Institute for providing data. R. Aanes was partly funded by the Norwegian Research Council through the KlimaEffekt programme (project no. 155903/720).
References Cited
FIGURE 1.
Brøgger peninsula, showing Ny-Ålesund, Austre Brøggerbreen (Br) Midre Lovénbreen (Ml), and natural barriers that impede migration from the peninsula (black lines). Inset shows location within Svalbard
FIGURE 2.
Population size N and log population growth rate G of Svalbard reindeer population on Brøgger peninsula, 1978–2002
FIGURE 3.
Winter mass balance for Midre Lovénbreen and Austre Brøggerbreen (1968–2002) averaged over the entire glacier (upper figure) and for the lowermost elevation intervals 50–150 m a.s.l. (lower figure)
FIGURE 5.
Regression coefficient r2 surfaces as a function of TSR and DDS, organized in rows by different values of TDDO, and in columns by results of R1, R2, and R3. Multiple surfaces in plot correspond to different values of DDI
FIGURE 6.
Histograms of four model parameters TDDO, TSR, DDS, and DDI (columns) obtained from the 95th percentile or better results of R1-R3 (rows)
FIGURE 7.
Histogram of number of rain or snow occurrences as a function of reported temperature at observation time, based on available data from Ny-Ålesund 1969–2002
FIGURE 8.
Raster plot of simulated winter snowpack (upper plots) and ground-ice thicknesses lower plots). Plots on left are for best-guess SBG and IBG (left panels) and on the right for the average of the 300 best simulations
FIGURE 9.
Winter mean ground-ice thickness ĪwBG compared to G (black points). Thick and dashed lines are parabolic and linear fits of ĪwBG to G. Gray points and lines show means and 1 S.D. of winter mean ground-ice thicknesses obtained from 300 best simulations
FIGURE 10.
May 1 snowpack thickness SM1BG compared to observed mass balance (black points). Line shows linear fit passing through the origin. Gray points and lines show means and 1 S.D. of May 1 snowpack thicknesses from 300 best simulations
FIGURE 11.
G observed (black line) and calculated using two regressions (3 and 6 in Table 3a). One regression (crosses) uses just [ĪwBG]2 and the second (circles) includes SM1BG and the previous year's population N(t − 1). The lower panel shows the regression residuals
TABLE 1
All recorded observations of ground-ice formation in Ny-Ålesund, and values assumed for use in R3
TABLE 2
Range and step-size (see text) of parameter values used in modeling
TABLE 3
Results of regressions of observed G against various parameters (see text) for the period 1978/79 to 2001/02. All regressions include a constant term. Table (a) uses all data, (b) removes data from 1993/94, and (c) also removes data from 1995/96. Each regression is numbered and includes one or more parameters (indicated by crosses). Bold indicates regression results significant at P = 0.05 or better
